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The equilibrium thermodynamics of the two dimensional Su-Schrieffer-Heeger Model is derived 
£NJ . by means of a path integral method which accounts for the variable range of the electronic hopping 

» I ' processes. While the lattice degrees of freedom are classical functions of time and are integrated out 

f> |_ exactly, the electron particle paths are treated quantum mechanically. The free energy of the system 

■ and its temperature derivatives are computed by summing at any T over the ensemble of relevant 

particle paths which mainly contribute to the total partition function. In the low T regime, the 
heat capacity over T ratio shows an upturn peculiar of a glassy like behavior. This feature is more 
sizeable in the square lattice than in the linear chain as the overall hopping potential contribution 
l ' to the total action is larger in higher dimensionality. 
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I. INTRODUCTION 

> \ 

, The lattice dimensionality is a key parameter in ma- 
& ■ terials science. One dimensional (ID) systems with half 
£h I filled band undergo a structural distortion [1] which in- 
i . creases the elastic energy and opens a gap at the Fermi 
surface thus lowering the electronic energy. The competi- 
tion between lattice and electronic subsystems stabilizes 
£j ' the ID structure which accordingly acquires semicon- 
'— — 1 , ducting properties whereas the behavior of the 3D system 
would be metallic like. Conjugated polymers, take poly- 
acetylene as prototype, show anisotropic electrical and 
optical properties [2] due to intrinsic derealization of 
I/"") ' 7r electrons along the chain of CH units. As the intra- 
CNJ , chain bonding between adjacent CH monomers is much 
stronger than the interchain coupling the lattice is quasi- 
1D. Hence, as a result of the Peierls instability, polyacety- 
lene shows an alternation of short and long neighbor- 
ing carbon bonds, a dimerization, accompanied by a two 
fold degenerate ground state energy. The Su-Schrieffer- 
| Heeger (SSH) model Hamiltonian [3] has become a suc- 
cessful tool in polymer physics as it hosts the peculiar 
'"O ' ground state excitations of the ID conjugated structure 
and it accounts for a broad range of polymer properties 
O ■ [4]. As a fundamental feature of the SSH Hamiltonian 
the electronic hopping integral linearly depends on the 
relative displacement between adjacent atomic sites thus 
leading to a nonlocal e-ph coupling with vertex function 
depending both on the electronic and the phononic wave 
vector. The latter property induces, in the Matsubara 
formalism [5] , an electron hopping associated with a time 
dependent lattice displacement. As a consequence time 
retarded electron-phonon interactions arise in the system 
yielding a source current which depends both on time 
and on the electron path coordinates. This causes large 
e-ph anharmonicities in the equilibrium thermodynam- 
ics of the SSH model [6]. Hopping of electrons from site 
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to site accompanied by a coupling to the lattice vibra- 
tion modes is a fundamental process [7] determining the 
transport [8] and equilibrium properties [9] of many body 
systems. A variable range hopping may introduce some 
degree of disorder thus affecting the charge mobility [10] 
and the thermodynamic functions. 

This paper focusses on this issue, dealing with the ther- 
modynamical properties of the SSH model in two dimen- 
sions and comparing them with the results obtained in 
one dimension [11]. Only a few extensions of the SSH 
Hamiltonian to higher dimensionality appear in the lit- 
erature [12] mainly concerning the phase diagrams [13] 
and the ground state excitations [14-16]. We apply a 
path integral method [17] which fully accounts for the 
time retarded e-ph interactions and, exploiting the above 
mentioned Hamiltonian linear dependence on the atomic 
displacement, allows us to derive the electron-phonon 
source action in two dimensions. The general formalism 
is outlined in Section II while the results are reported on 
in Section III. The conclusions arc drawn in Section IV. 



II. PATH INTEGRAL FORMALISM 

In a square lattice with isotropic nearest neighbors 
hopping integral J, the SSH Hamiltonian for electrons 
plus e-ph interactions reads: 

H= YX^Mft+lJr.s+h.c) 
r.s 

+ (Jr,s)y{fl, s + ifr,s + /l.C.) 

{J r ,s)x= -^[J-aAu x ] 

(J r ,s)y = ~[j-aAUy] 

Au x = u x (r + 1, s) - u x (r, s) 
Au y — u y (r, s + 1) - u y (r, s) 



(0.1) 



where a is the electron-phonon coupling, u(r, s) is the 
dimerization coordinate indicating the displacement of 
the monomer group on the (r, s)— lattice site, f£ s and 
f r _ s create and destroy electrons (i.e., 7r band electrons 
in polyacetylene) . The phonon Hamiltonian is given by a 
set of 2D classical harmonic oscillators. The two addenda 
in (1) deal with one dimensional e-ph couplings along 
the x and y axis respectively, with first neighbors elec- 
tron hopping. Second neighbors hopping processes (with 
overlap integral J^) may be accounted for by adding to 
the Hamiltonian the term such that 

# (2) - (Jr,-)x,,,(/ r t + l,.+l/r,- + '».C.) 

(Jr,s)x, v = ~\ J (2) - a^(Au x ) 2 + (Auy) 2 

(0.2) 
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The real space Hamiltonian in (1) can be transformed 
into a time dependent Hamiltonian [18] by introducing 
the electron coordinates: i) (x(t'), y{r')) at the (r, s) lat- 
tice site, ii) (x(t), u(t')) at the (r + l,s) lattice site and 
iii) (x(t'), y(r)) at the (r, s + 1) lattice site, respectively, 
r and t' vary on the scale of the inverse temperature (3. 
The spatial e-ph correlations contained in (1) are mapped 
onto the time axis by changing: u x r y \(r, s) —* u x ( y \{r'), 
u x {r+\, s) — > u x {t) and u y (r, s+1) — » u v (t). Now we set 
r' = 0, (^(r'),y(r')) = (0,0), (u x (t'), %(t')) = (0,0). 
Accordingly, (1) transforms into the time dependent 
Hamiltonian: 

H(t) = J x (r)(/t( a; (T),0)/(0,0) + / i .c.) 

+ Jy(r)(.fHo,y(r))f (0,0) + h.c) 

Jx{t) = -\[ J - <xu x (t)] 

J v {t) = -i[j-OUy(T)] 

(0.3) 

While the ground state of the ID SSH Hamiltonian 
is twofold degenerate, the degree of phase degeneracy is 
believed to be much higher in 2D [19] as many lattice dis- 
tortion modes contribute to open the gap at the Fermi 
surface. Nonetheless, as in ID, these phases are con- 
nected by localized and nonlinear excitations, the soliton 
solutions. Thus, also in 2D both electron hopping be- 
tween solitons [20] and thermal excitation of electrons to 
band states may take place within the model. These fea- 
tures are accounted for by the time dependent version of 
the Hamiltonian. 

As t varies continuously on the (3 scale and the r- 
dependent displacement fields are continuous variables 
(whose amplitudes are in principle unbound in the path 
integral), long range hopping processes are automatically 
included in H (r) which is therefore more general than the 
real space SSH Hamiltonian in (1) (and (2)). Thus by 
means of the path integral formalism we look at the low 
temperature thermodynamical behavior both in ID and 
2D searching for those features which may be ascribable 
to some local disorder related to the variable range of the 
hopping processes. 

The semiclassical nature of the model is evident from 
(3) in which quantum mechanical degrees of freedom in- 
teract with the classical variables u x i y \(T). Averaging 
the electron operators over the ground state we obtain 
the time dependent semiclassical energy per lattice site 
N: 

< H(t) > 



N 



J x (t)P(J, t, x(t)) + J v (t)P(J, t, y(r)) 

P(J, t, v(t)) = ^ J dkcos[k • v(t)] cosh(e k r)n F (e k ) 

(0.4) 

with v(t) = (x(t),0) and v(t) = (0,y(r)) in 
the first and second addendum respectively. — 
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—J^2 i=x y cos(fci) is the electron dispersion relation and 
up is the Fermi function. Eq.(4) can be rewritten in a 
way suitable to the path integral approach by defining 



< H{t) > 



N 

V(x(r)) 
V(y(r)) 
j(v(r)) = 
u(r) = (u x (t),u v (t)) 



= V(x(r))+V(y(T))+u(T)-i(v(r)) 

-JP(J,t,x(t)) 
-JP(j,r,y(T)) 
aP(j,r,v(r)) 



(0.5) 

V(x(t)) and V(y(r)) arc the effective terms account- 
ing for the t dependent electronic hopping while j(v(r)) 
is interpreted as the external source [21] current for the 
oscillator field u(t). Averaging the electrons over the 
ground state we neglect the fermion-fermion correlations 
[22]. This approximation however is not expected to af- 
fect substantially the following calculations. Taking a 
bath of N 2D oscillators, we generally write the SSH 
electron path integral, ((t) = (x(r),j/(r)), as: 

N 



<C03) I C(0) >= ft / Dudr) f D((t) 

]^r^(u, 2 (r)|^(r) 



exp 



exp 



N 



dr[^C 2 (r) + V(x(T))+V(y(r)) 



+ 5>(t).j(v(t)) 



(0.6) 



where m is the electron mass, M is the atomic mass 
and Ui is the oscillator frequency. As a main feature we 
notice that the interacting energy is linear in the atomic 
displacement field. Then, the electronic path integral 
can be derived after integrating out the oscillator degrees 
of freedom which are decoupled along the x and y axis. 
Thus, we get: 



N 



<c(/?)ic(o) >=n^ 



/ 



Dx{t) 



■ exp 



dr[-x 2 (r) + V(x(r)))--A(x(r)) 



A(x(r)) 



N 



i 



drj(x(T)) 



4M ^ ftwi smh(hwi/3/2) 



dr" cosh\Jiuji(\T 

1 2 



(3/2))j{x{r")), 



2sinh(nWi/3/2) 
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(0.7) 



and the 2D electron path integral is obtained after 
squaring the sum over one dimensional electron paths. 
This permits to reduce the computational problem which 
is nonetheless highly time consuming particularly in the 
low temperature limit. Note infact that the source cur- 
rent j(x(r)) requires integration over the 2D Brillouin 
Zone (BZ) according to (4) and (5) and this occurs 
for any choice of the electron path coordinates. The 
quadratic (in the coupling a) source action A(x(t)) is 
time retarded as the particle moving through the lattice 
drags the excitations of the oscillator fields which take a 
time to readjust to the electron motion. When the in- 
teraction is sufficiently strong the conditions for polaron 
formation [23] may be fulfilled in the system according 
to the degree of adiabaticity [24]. However the present 
path integral description is valid independently of the ex- 
istence of polarons as it applies also to the weak coupling 
regime. Assuming periodic conditions x(t) = x(t + (3), 
the particle paths can be expanded in Fourier compo- 
nents 

oo 

x(t) = x Q + 2^3?x„ cos((jJ„t) — 3x„ sin(w n T)^ 

71=1 

u n = 2im//3 (0.8) 

and the open ends integral over the paths J Dx(t) 
transforms into the measure of integration § Dx(t). Tak- 
ing: 

(0.9) 

which normalizes the kinetic energy term in (7), we 
proceed to integrate (7) in order to derive the full parti- 
tion function of the system versus temperature. 



f°° d$tx„ f°° d^sx 

J — oo J — oo 

(ttH 2 K bT j ma;?,) 



III. THERMODYNAMICS IN ID AND 2D 

Taking a dimensionless path (in units of the lattice 
constant a — IA) with: a n = 23?x„ and b n = — 2Qx n the 
functional measure in (9) can be rewritten as: 

t Dx(r) « (2Np+1) n( 2 ™) 2 / d -o 

J (2A m J „=i J -a 

2A f2A 

da n / db n 

-2A J-2A 



where N p is the cutoff on the Fourier components in 
(8). At a given temperature and for a given set of particle 
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path coefficients we first carry out the momentum inte- 
grations required by (4) summing over 1600 points in the 
reduced 2D BZ. This prepares the source current j(x(r)) 
to be used in the double time integration which yields, 
according to (7), the source action A(x(t)). Afterwards 
we sum the exponential of the total action over an ap- 
propriate set of particle paths. Numerically stable results 
are achieved by tuning the number of pairs N p in (10), 
the cutoff (A) on the integration range of the Fourier co- 
efficients in (10) and the related number of points (N\) 
in the measure of integration. N p = 2 suffices in the 
Fourier expansion of the path. The argument to set the 
cutoff A on the integration range is based on the fact 
that the functional measure normalizes the kinetic term 



From (11) we get, A <~ 3A m /\/27r 3 hence, A scales ver- 

sus temperature as A cx 1/VT. Computing (7) with (10) 
we find that infact a smaller cutoff, A ~ A m / (10\/27r 3 ), 
guarantees convergence in the interacting partition func- 
tion. Accordingly the number of points in each inte- 
gration range has to increase by decreasing T as N\ cx 
1/VT. In the 2D problem, A A ~ 35/ VT points for each 
Fourier coefficient are required in (10). Then we carry out 
the total 2D path integral by evaluating the contribution 
of (N\ + l) 2jv p+ 1 paths selecting the integer part of N\ 
at any temperature. Computation of the second order 
derivatives of the free energy in the range T G [1, 301]AT, 
with a spacing of 3K, takes 55hours and 15 minutes on a 
Pentium 4. 

The physical input parameters characterizing both the 
ID and 2D models are: the bare hopping integral J, 
the oscillator frequencies Wj and the effective coupling 
X = a 2 h 2 /M (in units meV 3 ). We take here a narrow 
band system with J = lOOmeV. The crossover between 
weak and strong e-ph coupling is given in our model by 
Xc ~ 7rJ7i 2 u4 ax /64 [2]. where uj max is the largestenergy 
in the oscillators bath. Let's assume a bath of N = 10 
low energy oscillators with huji <~ 20meV thus setting 
the crossover value at Xc ~ 2000meV 3 . A higher num- 
ber of oscillators in the same energy range would not sig- 
nificantly modify the results hereafter presented whereas 
lower oji values would yield a larger contribution to the 
phonon partition function mainly at low T. In Fig.l, a 
comparison between the ID and the 2D free energies is 
presented for two values of Xi one lying in the weak and 
one in the strong e-ph coupling regime. The oscillator 
free energies F p h are plotted separately while the free en- 
ergies arising from the total action in (7) , shortly termed 
Fsource, are due to the purely electronic plus e-ph con- 
tributions. In general, the 2D free energies have a larger 
gradient (versus temperature) than the corresponding ID 
terms. The F p h lie above the F source both in ID and 2D 
because of the choice of the hivf, lowering the energy of 



in (7): 




(0.11) 
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some oscillator modes would produce a crossing point 
between F p h and F source with temperature location de- 
pending on the value of \. Note that the F source plots 
have a positive temperature derivative in the low tem- 
perature regime and this feature is more pronounced in 
2D. Infact, at low T, the source action is dominated by 
the hopping potential V(x(t)) while, at increasing T, the 
e-ph effects become progressively more important (as the 
bifurcation between the Xi and xi curves shows) and the 
Fsource get a negative derivative. In 2D, the weight of the 
V{x(t)) term is larger because there is a higher hopping 
probability. This physical expectation is taken into ac- 
count by the path integral method. At any temperature, 
we monitor the ensamble of relevant particle paths over 
which the hopping potential is evaluated. For a selected 
set of Fourier components in (10) the hopping decreases 
by lowering T but its value is still significant at low T. 
Considered that: i) the total action is obtained after a dr 
integration of V(x(t)) and ii) the dr integration range is 
larger at lower temperatures, we explain why the overall 
hopping potential contribution to the total action is re- 
sponsible for the anomalous free energy behavior at low 
T. 

In Fig. 2, the heat capacity contributions due to the 
oscillators (C p h) and electrons plus e-ph coupling (C sou ) 
are reported on. The values are normalized over N. The 
dimensionality effects are seen to be large and, for a given 
dimensionality, the role of the e-ph interactions is mag- 
nified at increasing T. The total heat capacity (C p h + 
Csou) over T ratios are plotted in Fig. 3 in the low T 
regime to emphasize the presence of an anomalous up- 
turn which appears at T X 10K in ID and T < 20K in 
2D. This feature in the heat capacity linear coefficient is 
ultimately related to the above discussed overall hopping 
potential effects while the strength of the coupling \ has 
a minor role in the low T limit. 

Finally we emphasize that the degeneracy of the 
ground state scarcely affects the finite temperature equi- 
librium properties. This can be shown in ID by treating 
the SSH Hamiltonian in mean field approximation with 
a staggered dimerization amplitude which yields the size 
of the lattice distortion. By diagonalizing the Hamilto- 
nian in momentum space we see that the competition 
between electron and lattice energy determines the con- 
tribution to the partition function and shifts downwards 
the free energy versus temperature. Hence, the heat ca- 
pacity due to the degeneracy is very small (in compari- 
son with the values given in Fig. 2) and strictly linear in 
T. Accordingly, the peculiar and sizeable upturn in the 
heat capacity over T ratio does not depend on the nature 
of the ground state. The inclusion of the quantum lat- 
tice fluctuations beyond the mean field approximation is 
expected to reinforce our conclusion as the fluctuations 
tend to suppress the long range order associated with the 
dimcrized structure [25]. 
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IV. CONCLUSIONS 



We have investigated the thermodynamics of a semi- 
classical Su-Schricffer-Hecger (SSH) Hamiltonian ex- 
tended to two dimensions in which the lattice displace- 
ments are taken as classical variables. By mapping the 
real space interactions onto the time (r) scale we deal 
with an Hamiltonian accounting for a variable range of 
the electronic transfer integral. Such a model can be 
suitably treated in the path integral formalism as both 
the electron and lattice degrees of freedom are continuous 
functions of r which varies up to the inverse temperature. 
Although we are thus assuming a continuum model the 
present path integral calculation is valid for any range of 
e-ph coupling and for any value of physical parameters. 

The two dimensional path integral for the SSH par- 
ticle in a bath of harmonic oscillators has been derived 
by summing over a set of electron paths that becomes 
increasingly large at decreasing temperatures. Thus, the 
computation of the total action and of the free energy 
derivatives is highly time consuming in the low T regime 
where the quantum mechanical nature of the electronic 
subsystem is fully accounted for. At low T, the electron 
hopping integral largely shapes the free energy which dis- 
plays an anomalous positive derivative more pronounced 
in 2D than in ID. This behavior is mirrored in the upturn 
of the heat capacity over T ratio. For a given set of input 
parameters characterizing electron band, oscillator ener- 
gies and e-ph couplings, we find that the upturn shows up 
at T < 20K in the square lattice and at T ^ 10K in the 
linear chain. Also the size of the upturn is larger in two 
dimensions consistently with the fact that the 2D hop- 
ping potential contribution to the total action is higher 
than in ID. The presented results have been obtained as- 
suming a narrow electron band (J <~ lOOmeV) but the 
upturn in the heat capacity over T ratio appears also in 
the case of wide band systems. 

The degenerate ground state of the SSH model is not 
responsible for the anomalous behavior of the heat ca- 
pacity over T ratio. Infact, closing the linear chain in a 
ring geometry, one can easily diagonalizc the ID Hamil- 
tonian in mean field approximation and show that the 
contribution to the total heat capacity, due to the two 
fold degeneracy of the ground state, is small and linear 
versus temperature. A quantum description of the lat- 
tice displacement fields should further reduce the role of 
the dimerized ground state as the fluctuations tend to 
suppress the order parameter. Similar conclusions are 
expected to hold in the 2D case where, however, the de- 
generacy of the ground state is manyfold. 

In conclusion, while in the low T regime the thermo- 
dynamics is essentially governed by the electronic term 
of the Hamiltonian, at increasing T the strength of e-ph 
coupling progressively determines the slope of the free 
energy and heat capacity. Our results suggest that a 
glassy like thermodynamical behavior is a peculiar fea- 
ture of Hamiltonian models with a variable range for the 
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electronic hopping. Such a behavior is more likely to 
occur in higher dimensionality as the structure is more 
closely packed and, accordingly, the hopping probability 
is larger. 



FIG. 1. Colour online. Phonon {F p h) and Source Term 
(F sou ) contributions to the ID and 2D free energies for two 
values of the effective coupling x- Xi = 1440meV 3 (weak 
e-ph coupling) \2 = 2560me V 3 (strong e-ph coupling) . 



FIG. 2. Colour online. Phonon and Source Term contribu- 
tions (normalized over the number of oscillators) to the ID 
and 2D heat capacities for the same parameters as in Fig.l. 
The oscillator heat capacities are also plotted. 
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FIG. 3. Colour online. Total heat capacity over tempera- 
ture for the same parameters as in Fig.l. 
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